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Abstract 

We examine the conformations of a model for a short segment of closed 
DNA. The molecule is represented as a cylindrically symmetric elastic rod 
with a constraint corresponding to a specification of the linking number. We 
obtain analytic expressions leading to the spatial configuration of a family of 
solutions representing distortions that interpolate between the circular form 
of DNA and a figure-eight form that represents the onset of interwinding. We 
are also able to generate knotted loops. We suggest ways to use our approach 
to produce other configurations relevant to studies of DNA structure. The 
stability of the distorted configurations is assessed, along with the effects of 
fluctuations on the free energy of the various configurations. 
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I. INTRODUCTION 



If a segment of untwisted rod is forced to close smoothly on itself, it will take on a 
circular shape in order to minimize the elastic energy. Thermal fluctuations will exert a 
significant effect on this configuration only if the circumference of the circle is greater than 
the temperature-dependent persistence length of this segment. When the looped rod is also 
forced to undergo a twist, the circular shape may prove to be unstable. The rod will then 
either distort into a non-planar, non-intersecting form, or it may wind about itself, in an 
interwound, or "plectonemic," configuration PH7J. 

When torsional and flexing stresses are not too severe, || the DNA molecule can be 
modeled as such a rod. A short enough segment of looped DNA will thus take on shapes 
determined by the minization of elastic energy. Planar shapes, characteristic of "relaxed" 
DNA have been observed, as have the plectonemic forms when the molecule is over-or 
underwound. The distortion, or supercoiling, of DNA under the influence of torsional stresses 
is widely believed to have significant implications with respect to the action of this molecule 
in biological settings . That this is so may be inferred from the existence of enzymes known 
as topoisomerases, which act to facilitate the alteration of the torsional characteristics of 
DNA configurations 

The elastic model of DNA has been the subject of intense research in the past 30 years. 
Theoretical appraches to this problem include Lagrangian mechanics [jl3| — pTT|l , b-splines and 
(numerical) molecular dynamics ||18|| , and statistical mechanics [flQ ] . Only recently has it 
been possible to determine the equilibrium configuration of a closed segment of DNA, as 
modeled by distorted rod . A comprehensive discussion of the properties of these equi- 
librium configurations has not yet been published. 

In this paper we contribute to the knowledge of this configuration by describing a method 
by which one can obtain the shape of a closed section of an elastic rod, along with various 
key properties of this system, including its twist, writhe, linking number and energy. We 
will use this formalism to generate both unknotted and knotted forms of closed DNA. We 
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will also outline how our formalism can be exploited to produce other shapes of interest to 
DNA researchers. 

The configurations that will be discussed include those that arise as a result of the 
the supercoiling instability for a closed loop. We will also examine the configurations that 
represent the evolution of higher energy deformation of a circular loop. In addition, we will 
briefly discuss configurations that are associated with knotted loops. Such knots include, but 
are not limited to, the trefoil. All the more elaborate configurations are unstable with respect 
to fluctuations about the solutions that extremize the energy of the loop. However, it is 
possible to envision mechanisms, analogous to those leading to the nuclesome structure p|.|2T 



that stabilize these configurations and, therefore, cause them to be relevant to biological 
systems. 

As noted about, we will discuss the mechanical stability of the configurations to be 
studied, and will present results with regard to the contribution to the free energy of thermal 
fluctuations about the minimizing shape. These results will be relevant as long as thermal 
fluctuations do not give rise to a significantly alteration of the shape that the rod takes, and, 
as noted above, if the circumference of the rod is not large compared to the rod's persistence 
length. The precise criterion for this will be examined later on in this paper. 

The most important advance in the work presented here is the algorithm for the gen- 
eration of various closed structures and the analysis of the stability of the configurations, 
particularly the one associated with the principal supercoiling instability. On the basis of 
the stability analysis, it is clear that only some of the configurations that one can generate 
with the use of the methods described in this paper are, in and of themselves, mechan- 
ically stable. However, the existence of devices that support non-interwound supercoiled 
configurations in the case of open strands of DNA allow one to speculate as to the 

possible relevance of those configurations to DNA in biological settings. Furthermore, one 
can envision utilizing these methods in conjunction with other calculational tools to study 
the mechanical properties and packing configurations of longer strands of DNA. 

The paper is organized as follows. First, the mechanical model is described and impor- 



tant parameters are defined. Of especial importance is an identification of the parameters 
that will be adjusted in order to obtain solutions of the energy extremum equations. Key 
quantities, including topological invariants, will be defined in terms of those parameters. The 
stability of the circular configuration of a closed rod under torsional stress is then assessed. 
The principal instability is identified, and it is shown how to obtain the family of closed 
curves that represent the evolution of this instability. This family interpolates between the 
circle and a figure eight, that latter of which lies at the onset of interwinding. The stabil- 
ity of the members of this class is then assessed. It is found that they are, indeed, stable 
against small fluctuations. Then, it is shown that the method by which the outcome of the 
supercoiling instability is investigated can be utilized to generate deformed configurations, 
among which are configurations in which the closed loop is also knotted. Finally, there is 
a brief discussion of the way in which the contributions of fluctuations about the classical 
solution configurations. Details and background are relegated to appendices. 



II. ELASTIC MODEL AND EXPRESSIONS OF INTEREST 



In previous work we have outlined some developments in the elastic model of DNA |22 . 
The molecule is represented as a slender cylindrical elastic rod. At each point s the rod is 
characterized by relating the local coordinate frame T to the frame Tq rigidly embedded in 
the curve in its relaxed configuration. The relationship between the stressed and unstressed 
local frames is specified by Euler angles 9(s),<p(s),ip(s) needed to rotate T§ into T. 

We begin by summarizing some of the results needed for this paper. We denote the 
elastic constants of bending and torsional stiffness by A and C, respectively. The length of 
the rod is denoted by L. The elastic energy is given by 

E el = £ds (<j) 2 sin 2 6 + 6 2 ) +^(0cos# + ^) 2 ) (2.1) 

The "dot" notation for derivatives is used interchangeably with the explicit derivative with 
respect to arc-length, d/ds. 
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The twist is given by 



)- J ds (</> cos # + (2.2) 



Tw 

A relevant result from our previous work is that Ck, the linking number, can be written as 

Ck = J-1 + -L f L U + ^) ds (2.3) 
2ir Jo K ' 



Where we have used Fuller's theorem p3[ to obtain the writhe, and White's theorem |24 
for the total link. The quantity X is an integer that encodes the knotting of the loop. It is 
equal to the number of signed crossings of the loop, assuming that it has been flattened out 
as much as possible. 

III. CANDIDATES: THE ELASTIC EXTREMA 
A. Stability and Excitations of the Circular Loop 

As a precursor to the discussion of the possible deformations of a closed loop, we review 
the family of excitations of a closed, circular loop, and the conditions under which this 
configurations may become unstable to one or more of these excitations. Detailed discussion 
of these issues has been relegated to Appendices [A] and |B[ The Appendices also contain 
derivations of key formulas presented below. 

The excitation of a circular loop is describable in terms of the deviation of the Euler 
angle 9 from the constant value of it/ 2 that it has in the case of a circle. If we write 
6(s) = tt/2 + j(s), then small excitations about a circular shape will be of the form 

/27m \ 

7(s) oc cos I — - — ho) (3.1) 

where n is an integer (possibly 0) and 5 is an arbitrary phase factor. The energy of such a 
sinusoidal excitation proportional to 



The requirement that the fluctuations leave the loop closed and, also, that they not 
change the linking number of the loop rule out n = and n = 1 as possibilities. The allowed 
deformation with the smallest energy has n — 2. Substituting this value into Eq. ( |3.2| ), we 
find for the energy cost of this excitation 

*-i(?)>-(iH <-> 

It is clear that this deformation of the circle can give rise to a lowering of the energy This 
state of affairs holds if 

Ck > ^V3 (3.4) 

In fact, the condition in Eq. ( |3.4j ) is just the requirement for the supercoiling instability of 
the closed loop. Deformations associated with higher values of n lower the energy of the 
loop when the condition 

.4 



Ck > — y/n 2 - 1 (3.5) 
holds. 

It is clear that the instability with the earliest onset, and the one that will prove dynam- 
ically "strongest" is the one with threshold as given by 



B. A Family of Curves 

We start by looking for a family of closed solutions to the Euler-Lagrange equations that 
represent the evolution of the deformation associated with the supercoiling instability, i.e. 
n = 2 in Eqs ( |3.1|) and (|3.2| ) above. This family of curves can be indexed by the writhe, Wr, 
of a member, which ranges from Wr = for a simple circle to Wr = 1 of the other limiting 
member of the family, a "figure eight." The actual constraint imposed on the closed curves 
that are of relevance to this discussion is that the link is fixed at a predetermined value. We 
shall therefore examine the conditions under which the family of writhing curves can satisfy 
the imposition of a linking number 
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C. Euler-Lagrange Minimization 



Our goal is to extremize elastic energy and keep the curve closed. The functional to 
consider for closed configurations of DNA is 

H= f L dE el - Ftz ds (3.6) 
Jo 

The first term of ( |3.6| ) is clear — the elastic energy must be extremized. The second term 
enforces the constraint that the curve minimizing the elastic energy also closes upon itself. In 
general, one expects to introduce Lagrange multipliers controlling the extent of the curve in 
all directions. However, one can concentrate on the net displacement in only the z direction, 
because it is always possible to transform to a a reference frame in which the curve is closed 
in the XY plane. Written explicitly as a function of Euler angles ( |3.6|) becomes 

rL a 



n= f ds ±U^^e + e 2 ) + 

Jo 2 v ' 



C 2 

+— Ucos6 + i)) -FcosO (3.7) 



2 

The extrema are found by applying Euler-Lagrange equations to fl3.7| ). Denoting the con- 
served quantities as = dH/d(p and = dH/dip we obtain 



J<t> — Jip cos 9 



(3.8) 



^ = - <pcos6 (3.9) 

The equation for 9 is a quadrature obtained by integrating &Hjdd = dj ds&Hjdd with E 
as the constant of integration. Defining u = cos 9 we see that the behavior of solutions is 
governed by a cubic polynomial in u: 



2 - 2 { 3 ^ (Eo -Fu)-±(J} + J}- 2J^u) 



u 



A v u ' A 2 
IF 

= —{u-a){u-b)(u-c) (3.10) 



where we order the roots c < u(= cos 0(s)) < b < a. Eq. ( |3.1U| ) requires that cos^(s) 



oscillates between c and b. For the distortions studied here, there are two complete ocillations 



as the the closed curve is traversed once. All the relevant quantities, including the shape of 
the curve, can be obtained using ( |3.10| ). For example: 



«( s ) d(j)ds , 

du with 



lu(s=o) ds du 
ds f~A 



(3.11) 



V IF ^ u - a )(u- h) ( u - c) 
The key to the solution of the equation for the distorted closed loop is the determination of 
the parameters F, a, b and c. 

D. Constraints 

To complete the analysis we must determine the parameters generated by the minimiza- 
tion procedure. One choice of parameters is the set of invariants and J^, the constant 
E , and the Lagrange multiplier F. However, as indicated above, a better practical choice is 
the Lagrange multiplier F and the parameters a, h and c. To determine the parameters we 
impose constraints on the curve. Figure [l] is a good visual guide to the geometric meaning 
of the constraints. They are 

A The loop closes on itself in the ^-direction. Because the curve consists of four segments, 
in each of which the variable u(s) goes from b to c, or from c to b, we have 

,L/4 

/ u{s)ds = (3.12) 

J 

B The loop goes through one quarter of a turn in each segment. This reduces to 

/•V 4 d6(s) 7T , 

/ —i-^ds = - 3.13 
Jo ds 2 K ' 

C The linking number takes on a predetermined value for the loop. This requirement 
leads to the following mathematical constraint on the solution to the Euler Lagrange 
equations. 

' L *£@-ds = Ck (3.14) 
ds 



S 



D Finally, the length of the loop takes on a predetermined value. Mathematically, 

[ L ds = L (3.15) 
Jo 

The above condition is not as tautological as it seems. The actual parameterization of 
the curve will be in terms of the dependence of quantities on the cosine of the Euler 
angle 9(s). 

E. Reparameterization 

At this point we reparameterize the problem in terms of F, a, b, c instead of F, E , J^, J^. 
Parametrizing the problem by the roots of the polynomial is extremely advantageous: it 
makes the analytic manipulations more transparent; it also streamlines the computational 
tasks. The two sets of parameters are related in the following manner: 

E = F(a + b + c)) (3.16) 

[af [af 
■h = ]J—(pi±P2); J^ = \\ — (P1TP2) (3-17) 



with pm = [(c ± 1)(6 ± l)(a ± 1)] 



1/2 



The choice of branch (±) is imposed by the family of configurations sought. For circular 
DNA without intrinsic curvature, takes a — , and correspondingly a +fj Let us make 
some definitions that render the notation more transparent: 

Ab=(b-c) ; Aa={a-c) (3.18) 

(3,9, 



Employing (|3.11|) we rewrite the constraint equations, fl3.12|) , ( p,13| ), and (|3.15|) , in terms of 



the new parameters: 



1 the choice of a particular branch is a non-trivial procedure 
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aK (g) = AaE (q) 



(3.20) 



= -T- 
Aa 



(a + l)(6+l) 



(a-l)(6-l) 



1/2 



n 



A6 



-1-c 



c-1 



1/2 



TT ( ^ 

n ,g 



32A 



K(g) 



(3.21) 
(3.22) 



AaL 2 

Where i^, and II are complete elliptic integrals of the first, second and third kind, respec 
tively. The issue of the constraint on linking number will be left for later discussion. 



F. Solving for a, b, c and F 

We find that the optimal procedure for calculation of parameters appropriate to a solution 
is to utilize Q3.22j ) to eliminate F, then ( |3.20| ) and (|3.21|) eliminate Aa and c. Because the 
writhe, Wr, is a monotonically increasing function of Ab, we make use of this property to 
classify distinguish between the members of a family of curves. Once the constants F, a, b, c 
are determined, the desired solutions and all the relevant quantities are computed via elliptic 
integrals. For example, the explicit expression for 9(s) in the first quarter of oscillation is 
(we have inverted ( p.ll| )): 

4> and ip are obtained similarly from (|3.8|) and ( |3.9|) . 

Figure |2| displays the family of curves computed in the manner discussed above. Since the 
constraint equations involve elliptic integrals, finding a solution on a computer is virtu- 
ally instanteneous; analytically and computationally elliptic integrals are equivalent to, say, 
arcsin. 



G. Bounding Members: Circle and Figure Eight 



Let us check whether the initial member of our family, a circle with Wr = joins 
smoothly with the previously known stable family of twisted circles ||25|| . The circle corre- 
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sponds to Wr = 0. A circle in the XY plane the curve must have bo = cq = 0. ( |3.21|) now 
states that 



1 



Veto + 1 - Vo.q-1 



which in turn gives 
2nA 1 



J, 



(y/ao + 1 + Va-o - l) = —7— 



Combining 



) and ( |3.9|) to obtain the twist, Tw, of the circle (|3.25| ) gives 



(3.24) 



(3.25) 



(3.26) 



It is also of interest to compute the Tw of the figure eight, which, like the circle, can be 
performed virtually by inspection. The figure lies in the YZ plane, which forces <fi to behave 
as follows: (refer to Figure [I] for visualization) 



(3.27) 



Combining ( |3.27| ) and ( |3.8| ) forces J<^ 8 = J^, 8 = which immediately sets Tw 8 = 0. The 
value of Ab$ is easily determined from the fact that c§ = — 1 (this can be seen from the 
curve itself) necessitating a = 1. Then ( |3.21| ) yields Afe 8 = 1.6522 . . .. 



IV. LINKING NUMBER AND THE PLECTONEMIC TRANSITION 

We have found a family of writhing solutions that are the extrema of elastic energy. The 
writhe of the curves covers the range < Wr < 1. 2 However, the physical constraint 
imposed on the molecule is Ck, the linking number. The explicit expressions for Tw, Ck 
and Wr are easily obtained: 

2tt A 
Tw = -=K(q)(p 1 +p 2 )- (4.1) 



2 That the figure eight has Wr = 1 just before crossing can be seen from the shape of the curve. 
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2n 

Wr = ^=[-K (q) (pi+pa) 
V Aa 



1 + c \ — 1 — c / 1 — c \ 1 — c 



(4.2) 



and, using White's theorem [24 



Ck = Tw + Wr (4.3) 

With reference to Figure [3|, we can now describe what happens to a deformed loop of 
DNA as the linking number is increased. The loop remains a circle until ACk = y/3A/C. 
After that there are three possibilities. If A/C < 0.5 then the writhing family supports a 
steady increase to the linking number limit Ck = 1 of the figure eight and the molecule 
folds continuously until self-crossing occurs. Further increase of Ck presumably results in 
a plectonemic configuration. If 0.5 < A/C < 1 the one can expect the curve to distort 
continuously until the writhe achieves some value intermediate between zero and one. There 
is, at this threshold value of the writhe, a transition, almost certainly to an interwound form. 
When A/C > 1, none of the members of the writhing family support the neccesary linking 
number, and as soon as ACk exceeds the supercoiling threshold value \/3A/C the twisted 
circle snaps into a plectoneme. An interesting fact is that this behavior is independent of 
the length of the molecule. That this ought to be so follows from the absence of an absolute 
length scale in the problem of the deformed loop. 

Greater insight into the behavior of the deforming loop is gained if one also investigates 
the way in which the energy depends on the various properties of the loop. The energy can 
be expressed in terms of the previously introduced quantities. A straightforward calculation 
yields for the energy of the loop 

E = F^a + b + c) + - A ( Pl +p 2 ) 2 --iXjds (4.4) 

Where the final integral is over arclength. Expressing the circumference of the loop in 
terms of F, a, b and c, we obtain the following result for the energy in terms of the total 
circumference of the loop, /, the bending modulus, A, and the parmeters a, b, and c: 
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E = 2 A {(„ + i + e) + I ( Pl + [i - l] } (/ d " ) (4.5) 

1 ^ 4 LL JJ \ Jc J (a - u){b - u){u - c) j 

The energy of the loops as a function of writhe for various values of the modulus ratio 
A/C is plotted in Figure |. A few facts about the energies of the deformed configuration 
can be established numerically (and we do not doubt that analytic demonstrations can also 
be constructed). First, the energy is a monotonically increasing function of linking number, 
if not of writhe. This is not immediately evident in Figure |], although it is indicated, in 
that the energy increases or decreases monotonically with writhe when the link does so. 
Furthermore, when there is a maximum in link as a function of writhe, there is also an 
energy maximum, and, as indicated in Figures |3| and [|, the maxima occur at the same value 
of the writhe. 

In addition, when there is a maximum in the linking number as the writhe increases from 
zero to one, so that there are two possible writhes associated with the same linking number 
in for a range of the latter quantity, the configuration with the lower writhe always has the 
lower energy. This can be established by a separate set of calculations. This implies that 
when, for instance, 0.5 < A/C < 1, the loop, in distorting from its planar form, will not 
be susceptible to a discontinuous distortion to a more highly writhed, but not interwound, 
configuration before it "snaps" to a plectoneme. In the absence of thermal fluctuations or 
external perturbations, that transition will to occur when the link has reached its maximum 
value function of writhe. 

For a more detailed discussion of the energetics of the deforming loop, along with a 
calcuation of the energy of the plectoneme, the reader is referred to |20| 



V. STABILITY OF THE DEFORMED LOOP 

Given the existence of the solutions for the deformed loop corresponding to the evolution 
of the supercoiling instability, it is desirable to determine whether these solutions are, them- 
selves, stable against small fluctuations. There is no a priori guarantee that this is so. We 
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will find, in fact, that the stability of the closed configuration is enforced by the relatively 
large number of constraints that must be satisfied by fluctuations about the extremizing 
solutions to the Euler-Lagrange equations. A number of technical details associated with 
the determination of the stability of a configuration are to be found in various appendices. 



A. Stability in the Absence of Constraints: the Translation Mode 

The energy of a fluctuation about a solution to the Euler-Lagrange equation is expressible 
in terms of solutions to the linear second order differential equation (|A9|) , equivalent to the 
Schrodinger equation of a particle in the potential ( |A10| ). In particular, any fluctuation, ((s) 
of the Euler angle 8(s) about the form it takes in an extremizing solution can be written in 
the form 

COO = (5.i) 



where \&z(s) is a solution to Eq. ( |A9| ) , with eigenvalue A/. Assuming that the \EVs are 
normalized, the energy of this fluctuation is given by 

J2 K ? X i (5-2) 
i 

The extremizing solution will be stable as long as there are no solutions to Eq. (|A9|) with 



negative eigenvalues. Now, consider the potential V(s), displayed in Figure ||, associated 
with one particular member of the family of solutions that we have been considering. 

A striking property of this potential is that it is always negative. From this one can infer 
that there are solutions of the Schrodinger-like equation with negative eigenvalues. In fact, 
the existence of such states is guaranteed by the existence of the translational mode — the 
"fluctuation" equal to the derivative with respect to arc length of the extremizing 6(s). As 
the system is invariant with respect to translations along the closed loop, the infinitesimal 
transformation of the extremizing solution 8 c i(s) — > 6 c i(s + 5s) = 9 c i{s) + 5sd8 c i(s)/ds has no 
effect on the energy of the configuration. The function $t(s) oc d9 c i(s)/ds has the requisite 
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periodicity, in that it is unchanged if s goes to s±nL. This implies the existence of a solution 
to the Schrodinger equation with zero eigenvalue. This solution is displayed in Figure |6|. In 
both Figure |5] and |6], the arclength parameter ranges from —L/2 to L/2. This will prove 
to be useful later on, when we take advantage of the reflection invariance of the portential 
V(s). 

The important thing to note is that the translational mode has nodes. Note, furthermore, 
that this mode is odd on reflection about s = 0. On the basis of a simple node-counting 
argument, one can readily establish that there will be another antisymmetric solution to 
the equation (|A9|) , having fewer nodes in the interval —L/2 < s < L/2, and, hence, a lower 
eigenvalue than the translational mode. This lower eigenvalue is necessarily negative. As it 
turns out, there are two symmetric solutions to (|A9|) having the requisite periodicity and 
negative A's. The extremizing solution is, thus, nominally unstable with respect to three 
kinds of fluctuations. 

B. Effects of Constraints 

Certain constraints apply to any fluctuation in a closed loop. In fact, there are five 
such constraints, listed in Appendix |A]. These constraints must be incorporated into any 
calculation of the stability of an extremizing solution. The general effect of such constraints 
on the stability calculation is outlined in Appendix [Cj The determination of the stability 
of a deformed loop amounts to a search for zeros of a determinant obtained by sandwiching 
a response kernel between five different functions, each associated with one of the five con- 
straints that must be satisfied by any small variation of the solution to the Euler Lagrange 
equations for the extremal loop configuration. The five functions are readily extracted from 
the integrals in Eqs. ( |A16| ) ( p\22| ) . With the use of symmetry arguments, one can verify that 
the five-by-five matrix constructed out of the response kernel, which is exhibited in Appendix 
[B], reduces to a four-by-for matrix and a one-by-one matrix. The latter matrix consists of 
the expectation value of the response kernel with respect to the function associated with 
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closure in the y direction. This function is contained in the integral in Eq. ( |A2C| ) . 

The matrices are straightforwardly constructed, although the calculation is somewhat 
tedious. We find that all the configurations associated with the supercoiling instability that 
interpolate between the circle and the figure eight are mechanically stable in that the for 
no negative value of the parameter A does the determinant pass through zero. An example 
of the calculation of this determinant for a particular deformation of the circular loop is 
displayed in Figures [7|, || and ||. 



VI. GENERALIZATIONS 

A. Higher Order Deformations 

The solution that interpolates between the supercoiling instability and the figure eight 
deformation bordering on the regime in which the loop interwinds does not represent the 
only possibility for deformation. One can also generate deformed loops associated with the 
evolution of higher energy excitations of the circular loop. This is accomplished by altering 
the requirement on the change in the Euler angle over a complete "period" of the oscillation 
of the variable u(s) = cos 6(s), as it cycles between its limiting values of b and c. If, instead 
of asking that <p advance by ir/2 in half a period, one requires a change of n/n, then it is 
possible to generate a family of solutions to the Euler-Lagrange equations associated with 



the excitations having n periods around the circumference of the circular loop. Figure [LO 
displays such a distortion of the circular loop. Here, n = 5. The stability analysis of this 
this family of solutions is straightforward. Starting with the translational mode, one counts 
nodes and determines the number of solutions to the eigenvalue equation for fluctuations 
for which A must be negative. There are simply too many to allow for stabilization by the 
action of constraints. In the absence of external stabilizing mechanisms, such as histones or 
their equivalent, a loop will spontaneously distort out of this configuration. 
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B. Knotted Configurations 



It is possible to generate solutions to the Euler Lagrange equations for the extremization 
of the elastic energy that close and that are, in addition, knotted. Again, one simply alters 
the requirement on the way in which the Euler angle advances over the course of a "period" 
in the oscillation of the quantity u(s) = cos8(s) between its two limiting values of b and 
c. Imposing the requirement that the change in <ft(s) is equal to mn/n, where m and n are 
relatively prime, and m < n, one obtains knotted solutions. For example, setting n = 3 
and m = 2, one generates closed loops in the form of trefoils. An example of this solution 
to the Euler-Lagrange equations is displayed in Figure Alternatively, setting n = 5 
and m = 3, one obtains a knotted solution with fivefold symmetry. This solution is shown 



in Figure [12]. Both displayed loops are members of a family having the same topological 
characteristics. These closed extremal curves interpolate between nearly circular, "braided" 
loops and relatively strongly writhing forms, such as are displayed in the figures. 

It is possible to assess the stability of the solutions shown in the figure. In this case the 
constraints on fluctuations do not suffice to stabilize them against deformations that lower 
their elastic energy. Thus, in the absence of externally imposed stabilizing mechanisms, 
these extremal configurations are not mechanically stable. 



The trefoil solution has also been generated in finite element calculations [IS|, and knot- 
ted configurations of closed DNA have long been known to exist in vivo |M 



VII. ENTROPIC CORRECTIONS 

The analysis of the fluctuations about the classical solution described abve readily lends 
itself to a calculation of the entropic contributions to the partition function of a loop of DNA 
at finite temperature. Performing an expansion of the Euler angles about the form taken in 
a solution to the Euler Lagrange equation and retaining terms that are second order in the 
deviation of those angles about their classical values, one obtains the following expression 
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for the partition function of a closed DNA segment 

Z oc e -Wo,^o) J p 7 ( s ) exp [-/3< 7 (s)|£|7(s))] (7.1) 



where the operator £ is given by 



£=-^ + V(.) (7.2) 



with V(s) as defined in Eq. ( [A10| ). The quantity 7(5) is the deviation of the Euler angle 
6(s) from its classical value. 



The integral over 7(5) in Eq. ( |7. 1|) yields the inverse square root of the Fredholm deter- 
minant, F(0), of the operator C, where 

F(A)=n(^-A) (7.3) 

1 

The quantities A; are the eigenvalues of the operator C The Fredholm determinant is readily 
calculated with the use of a method commonly exploited in the study of instanton effects 
in nonlinear systems [32|. The application of this method to the case at hand is outlined 



in Appendix [EJ. One finds that the determinant can be expressed in terms of the quantity 
T(A), defined in Eq. (|E(f). This quantity is plotted for a characteristic member of family 
of deformed loops that interpolate between a circle and a figure eight in Figures [13] and [14]. 
A noteworthy property of the Fredholm determinant as displayed in the plots is the fact 
that this function of A possesses a zero at A = — this as a consequence of the existence 
of the translation mode — and the fact that there are zeros at negative values of A. These 
latter features point to the instability of the unconstrained solution to the Euler-Lagrange 
equations. 

Given that T(0), and, by extension, -F(O), is equal to zero, one expect the gaussian 
integration over the variable 7(5) to blow up. However, the translational mode has a well- 
defined influence on the fluctuation spectrum. It simply gives rise to a multiplicative factor 
reflecting the freedom one has in the fixing location the distortion on the loop. One eliminates 
this zero from the Fredholm determinant by noting that 



d 

dX 



ri(^-A) 



n (Ai - a,,,) 

^Z 



(7.4) 



This means that we can eliminate the translational mode from the Fredholm determinant 
by making the replacement 

dT{\) 



T(0) 



dX 



-T'(0) 



(7.5) 



A=0 



The result of an unconstrained gaussian integration over fluctuations is, then, proportional 
to 



p-W (-T'(0)r 1/2 



(7.6) 



where N is equal to the number of modes contributing to the fluctuation spectrum. In light 
of the behavior of T(X), as displayed in Figures O and [TJ], this result is clearly pathological. 



A. Constraints 

Of course, one is not allowed to integrate freely over all fluctuations. The constraints on 
the influence of those fluctuations are the same as applied in the linear stability analysis. 
The imposition of constraints leads to an additional multiplicative term 

(det(G(A = 0))r 1/2 (7.7) 

For the generaly principle underlying this result, see Appendix |F|. The matrix G is a 5 x 5 
version of the matrix introduced in Appendix y. The determinant of this matrix has poles 
at the negative values of A at which the function T(A) passes through zero. Figures [7], |8| 
and |9] display the dependence on A of the two quantities that when multiplied together yield 
this determinant. The plots in these figures are for the same solution to the Euler-Lagrange 
equation as generated the T(A) displayed in Figures [IB] and |TJ]. The quantity 

p-N/2 (_ T '( ) det (G (A = 0))r 1/2 (7.8) 

represents the contributions of fluctuations about the classical solution, to within readily 
calcuated numerical factors. In the example for which quantities are displayed, the argument 
of the square root is finite and positive. 
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VIII. CONCLUSION 



We have presented a formalism for obtaining the elastic mimima of a segment of closed 
DNA subject to a constraint in the linking number. The methods exploited here have been 
utilized to construct the family of deformations interpolating between the "relaxed" circular 
form taken by such a segment when the molecule is insufficiently under- or overwound to 
induce a supercoiling transition and the figure eight form that represents the threshold of 
interwinding. The members of this family are stable with respect to small fluctuations. The 
same methods also give rise to solutions that represent the evolution of "higher order" fluc- 
tuations about a circle of the twisted loop. These solutions to the Euler-Lagrange equations 
for the extremization of the elastic energy of a closed loop are unstable with respect to small 
fluctuations, in that there exist one or more fluctuation modes that lead to a lowering of 
the energy with respect to the classical solution. It is also possible to construct knotted 
solutions to the Euler-Lagrange equations. These configurations also represent saddle-point 
solutions to the extremization equation. However, as previously noted, it is possible to 
envision stabilizing mechanisms consistent with the known structure of biological sysems. 

In this paper, entire closed loops were generated and studied. The same method can, in 
principle also be utilized to construct a truly "finite" finite-element-analysis in which elastic 
models can be traced out with the use of non-infinitesimal segments. Given our experience 
in the project described herein, we are confident that the properties of the segments can be 
controlled and explicitly displayed. This ought to give rise to a substantial savings in time 
and effort in numerical studies of DNA structure. 
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APPENDIX A: STABILITY CONSIDERATIONS 



Starting with the expression for the elastic energy 




(Al) 
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where 



£ei(s) 



A 
~2 



'd6{s) 
ds 



+ sin 2 9(s) 



' d(f)(s) 
ds 



C 
~2 



dib(s) n/ s dd)(s) 
' + cos0(s) 



ds 



ds 



(A2) 



See Eq. ( |2.1| ). We expand the Euler angles about their extremizing values as follows 

<f>(s) = 4>d(s) + a(s) 

Hs) = Ms) + P(s) 
9(s) = 9 cl (s)+ 1 (s) 



(A3) 
(A4) 
(A5) 



The subscript cl stands for the extremizing, or "classical" values of the Euler angles. The 
deviations from these classical values, a(s), /3(s) and j(s), are assumed to be small. Substi- 
tuting from Eqs. (|A3|) - ( 5E ) into (|A2Q , one finds at quadratic order in the deviations from 
the extremizing solutions a local energy equal to 

At ■ r 

- Yi{s) 2 + d (s) 2 7 (s) 2 [cos 2 - sin 2 # c/ (s) 
+4d(s)7(s)0 c/ sin6> d (s) cos9 d (s) + d 2 sin 2 9 d (s 
$(s) + a(s) cos6 d (s) - <}> d (s)j(s) sm8 d (s) 
ipd(s) +<i> d {s) d(s)^(s)sm9 d (s) + 0(s)7 2 (s) cos0 d (s)/2 



C 
+ 2 



(A6) 



When it does not lead to confusion, the, the subscripts will be dropped from the "classical," 
or extremum, Euler angles. Because the classical Euler angles satisfy an extremum equation, 
there is no term linear in the deviations. 

A cursory investigation Eq. (|A6|) reveals that the second-order energy depends on a(s) 
and j3(s) only through their derivatives. Minimizing the energy with respect to these vari- 
ables, we are left with the following dependence of the "fluctuation" energy on the angular 
variable 7(5) 

A (d-f{s)\ 2 7 (s) 2 (4 + 4) (2cos 2 #(s) + 1) - cos 9 (cos 2 9 + 5) 



ds 



A sin 4 9 



(A7) 



Recall that the Euler angle 9(s) in Eq. ( |A7| ) is the solution to the minimzation equation. 
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The equation (|A7|) for the energy of a fluctuation can be further reduced if one makes 
use of the relationship between the quantities and and the roots a, b, and c of the 
cubic polynomial in Eq. ( |3.17| ). The new form of the energy is 



1 (*y(8)V 

2 V ds J 

1 (2 - u(s)) (1 - u(s)f (a + !)(&+ l)(c + 1) + (2 + u(s)) (1 + u(s)f (a - 1)(1 - 6)(1 - c) 2 
4 (l-u(s) 2 ) 2 1[S) 

(A8) 

The quantity u(s) satisfies Eq. (|3.10| ). This equation is equivalent to the expectation value 
of the energy of a particle in a one- dimensional potential. This expectation value can be 
expressed in terms of the eigensolutions and eigenvalues of the corresponding Schrodinger 
equation. The stationary version of the Schrodinger equation has the form 

_ d 2 V{s) 
ds 2 

where 



+ V(s)V(s) = Atf(s) (A9) 



V(s) = 

1 (2 - u{s)) (1 - u{s)f (o + l){b+ l)(c + 1) + (2 + u{s)) (1 + u{s)f {a - 1)(1 - 6)(1 - c) 

2 (1-m(s) 2 ) 2 

(A10) 

The stability of a fluctuation is tied to the sign of the eigenvalues of the stationary 
Schrodinger equation. If all the allowed values of A are positive, the solution about which 
fluctuations occur is positive. On the other hand if there is one or more negative A, then an 
instability exists. 

1. Constraints 

Fluctuations about the classical solution must obey certain constraints. In particular, 
they cannot change the following properties of the loop: 
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A. The loop closes smoothly: 

rL d(j)(s) 



o ds 



ds = 2nn (All) 



where n is an integer 

B. The net linking number is fixed: 

f L ^Ids = const. (A12) 
Jo ds 

C. The loop closes in the x direction: 

[ L sin 9(s) cos <f)(s)ds = (A13) 
Jo 

D. The loop closes in the y direction: 

[ L sin 9(s) sin (f)(s)ds = (A14) 
Jo 



E. The loop closes in the z direction: 

rL 



( cos6{s)ds = (A15) 
Jo 



In the above expressions, the Euler angles are not necessarily equal to their extremum values. 

Expanding the solution about its extremum form and expressing the fluctuations in the 
Euler angles and ip in terms of the fluctuation, 7(5), in the Euler angle 9(s), the conditions 
above take the forms shown below 



Jo n - n (s)) 3/2 v ' v 1 



A. Smooth closure: 

) +P2[U ' 

2AJo (i_ m ( s )) 3 / 2 

Here, and below, the quantity u(s) is cos^(s), where 9(s) is the solution to the ex- 
tremum equations. 
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B. Constancy of the link: 



F f M1 -u(s) f- M l + u(s) f ^ = 

Jo n - u (fiW 2 



2AJ (1-m(s) 2 ) 

C. Closure in the x direction: 

rL \ , x , ,x, fypi(l-u(s)) 2 + p 2 (l+u(s)) 2 



^ |w(s)cos0 d (s) + y — (i _ M ( s )2)3/2 2*(s) J 7(s) ds = (A18) 

Here, the quantity T x {s) is given by 

2*00 = f v /l- M ( S ') 2 sin0 c/ ( S , )rfs / (A19) 



D. Closure in the y direction: 



j o |w(s)sm0 d (s) - y — (i -^(5)2)3/2 V s ) I 7(s) = (A20) 



where 



l y (s)= I* Jl - u(s') 2 cos <p d (s') ds' (A21) 

JL/2 V 



E. Closure in the z direction: 



£ yf\ - m(s) 2 7 (s) rfs = (A22) 



APPENDIX B: STABILITY ANALYSIS OF THE CIRCULAR LOOP 



The solution of the extremum equations leading to a circular loop is 



u(s)=0 (Bl) 

Ms) = ^ (B2) 

27TS 

^ c/ (s)=£A;— (B3) 
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the quantity Ck being the linking number of the circular loop. Making use of Eq. 
find for the eigenvalue equation for fluctuations 



A 
~2 



d 2 l{s) 
ds 2 



Pcl{S) 

ds 



cy( dMs) ' 

A) [ ds 



7 (s) } = A 7 (s) 



, we 



(B4) 



substituting for the rate of change of the classical Euler angles, 4> c i(s) and (p c i{ s ), the eigen- 
value equation becomes 



A I d^(s) 
2 I ds 2 



2j^ 2 
L 



I : (-) Ck 2 



7 ( s ) = A 7 (s) 



(B5) 



Now, the solutions to the equation above are 

/27m 



"f{S) oc COS 



s + 5 



(B6) 



where n is an integer and 5 is an arbitrary phase angle. Looking at Eq. ( p5|) , one might be 
tempted to conclude that the circular loop is always unstable, in that solutions of the form 
( |B6|) with n = and n = 1 give rise to negative eigenvalues. However, those solutions are 
inconsistent with the constraints on fluctuations. A cursory inspection of the constraints 
listed in Eqs. ( |A16| ) to Q3.12j ) reveals that the following conditions must hold in the case of 
the circular loop 



o 



7 (s) ds = 



o 



7 (s) COS ( —r- ) <l* = 



V L 

/2ns 



~/{s sin . 
o V L 



ds = 



(B7) 
(B8) 
(B9) 



These three constraints on explicitly rule out flucutations of the form (|B6|) with n = or 
n — 1. All other values of n are allowed. If we substitute a fluctuation given by ( p6|) with 
n > 2, the equation for the eigenvalue A is 



A 



A /27TV 

2 \ l) 



n 2 - 1 



(BIO) 



The lowest value of A is associated with the smallest allowed value of n 2 , corresponding to 
n = 2. Replacing n by 2 in Eq. ( B10|) , we find 
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According to Eq. (pil| ), A will be negative, corresponding to an instability in the circular 
configuration, when Ck > y/3(A/C). 



APPENDIX C: GENERAL EFFECT OF CONSTRAINTS ON A STABILITY 

CALCULATION 

The question of the stability of a solution to the Euler-Lagrange equations is posed in 
terms of the eigenvalue spectrum of a linear operator. This, in turn, can be recast in terms 
of the problem of finding extremal values for the expectation value 

(£140 (Gi) 

where £ is the linear operator. The constraints are equivalent to requiring that the £ between 
which the operator is sandwiched is orthogonal to a set of m x's. There is also the constraint 
on the absolute magnitude of £. The constraints are, then of the form 

m = i (C2) 

(£lXi> = (C3) 

In Eq. flC3| ) , the index I runs from 1 to m. The equation for the extremum of the quadratic 
form (|ClD , subject to the constraints ( p2|) and (p3D , takes the form 



£|£> = Ai0+5>lx*> (C4) 
i=i 

The coefficients A and A/ are Lagrange multipliers, which enforce the constraints to which 
the system is subject. The solution to the above equation is 

m a 

10 = £7-41X1) (C5) 
1=1 L ~ A 

The Lagrange multipliers A; must now be adjusted to ensure the orthogonality requirements. 
These requirements are of the form 
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m j 

= J2MXkh r\xi) 

1=1 L ~ A 

= G^A, (C6) 

This set of m equations for the Lagrange multipliers A/ has non-trivial solutions only if the 
determinant of the m x m matrix G is zero. The equation \Gjk\ = represents a condition 
on the parameter A. 

Now, given a solution to Eq. ( |C6| ) , we take the expectation value Substituting 
from the right hand side of Eq. (|C6|) , we find for this expectation value 



171 A m A 771 A 

i=l L, — A l=1 L, — A l=1 L - A 

m 

= A (C7) 

In Eq. ( P7|) we have made use of the orthogonality of £ to the x/' s - We are also assuming 
that the function £ is normalized. Thus, in solving for the value of A that satisfies Eq. ( |C6[ ) 
we are also determining the effective values of the eigenvalues of the constrained problem. 



APPENDIX D: THE RESPONSE KERNEL 

The quantity 1/(C — A) represents the response kernel in the interval —L/2 < s < 
L/2. This kernel, which can be written in the form K(s,s'), is the inverse of the operator 
— d 2 /ds 2 + V(s) on that interval, and it has the additional property that it maps onto periodic 
functions. This response is constructed out of two solutions to the differential equation 

d 2 $(s) 

-^P + V(sMs) = ms) (Dl) 
The first solution, $i(s) is even under reflection about the origin. It has the property 

$x(0) = 1 (D2) 



ds 



= (D3) 

s=0 
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The second solution, $2(5), is odd with respect under reflection about the origin. It satisfies 
the following conditions 



$ 2 (0) = 



ds 



(D4) 
(D5) 



s=0 



That the response kernel below satisfies all requirements above can be established by 
explicit calculation: 

K(s,s') = $ 1 (s > )$ 2 (s < ) 



$ 2 (L/2) 



(D6) 



2$i(L/2) 

In Eq. (|D6|) , the argument s>(<) is the greater (lesser) of s, s', and the dots refer to 



differentiation with respect to arclength, s. 



APPENDIX E: THE CONSTRUCTION OF THE FREDHOLM DETERMINANT 

The Fredholm determinant of the operator £ — A, where 

£~£ + VM (ED 

is of the form 

II (Ai - A) (E2) 
1 

where the A^'s are the eigenvalues of the operator L. In the case of interest the eigenfunctions 
of this operator are periodic, in that if C^i = then \E^(s + L) — ^i(s). 

The Fredholm determinant is defined in terms of its analytic properties in the complex 
A plane. It is obvious that this function of the variable A has no poles at finite values of A, 
and that all zeros are on the real axis, at the locations of the eigenvalues of the operator C 
There is another function of A having this property, formed from the solutions <&i(s) and 

30 



$2(5), defined in Appendix [D[ To construct this function we note that any solution of the 
equation £^ = A^ can be represented in terms of the two linearly independent functions 
$1 and $2 as follows 



(s) = tf(0)$i(s) + 



ds 



(E3) 



s=0 



This means that we can express the function ^(s) and its derivative, ^(s) at s = L in terms 
of the function and its derivative at s = in the following form 



» ' (0) ^ 



V 



/ 



*(0) 



(E4) 



$ 2 (L) 

Given the fact that the Wronskian of $i and $2 is equal to one, the matrix on the right 
hand side of Eq. ( |Ejj ) has a determinant of one. Now, if the function ^(s) is periodic in s, 
then the right hand side of Eq. ( |E4| ) is equal to the right hand side. This leads immediately 
to the characteristic equation for the matrix 







$i(X) — 1 $ 2 (L) 

$ 2 (L) - 1 



= 2-$i(L) 
This means that the function of A 

T(A) = -2 + $x(L) + 



ds 



(E5) 



ds 



(E6) 



is equal to zero at the eigenvalues of the operator £ = —d 2 /ds 2 + V(s). Furthermore, it can 
be established that this function is free of singularities and zeros at any other finite value of 
A. 

The function T(A) is not identically equal to the Fredholm determinant of the operator 
£ — A, in that the behavior at |A| = 00 of the two quantities is not the same. However, if 
we define T (A) as the function corresponding to T(A) when the potential V(s) has been set 
equal to zero, and if we denote by F(X) and Fo(A) the corresponding Fredholm determinants, 
then the following relationship can be esablished 
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m_ = Tw (F7) 

Fo(X) T (A) 1 1 

Given that the Fredholm determinant of the operator £ — A = —d 2 /ds 2 — A is readily 
calculated, Eq. ( |E7| ) leads to a direct determination of the desired quantity. Figures [13] and 



14] display a characteristic T(A) for one of the writhing solutions that interpolate between a 
circle and a figure eight. Figure [14] shows in detail the behavior of the Fredholm determinant 
in the vicinity of A = 0. Note that this function passes through zero at A = 0, and that 
there is a zero of this function at a small, positive value of A. 

APPENDIX F: THE EFFECT OF CONSTRAINTS ON A GAUSSIAN INTEGRAL 



Given the quadratic form 



the gaussian integral 



XiAijXj = x ■ A ■ x (Fl) 



/n In \ 

Y[ dxi exp - x iAjXj (F2) 
t=i v ij=i J 



subject to m the constraints 

x-yi = 0, \<l<m (F3) 

is equal to 



m ~ m 



/n / m \ 

][J c?Xj exp —x • A • x + £ ik>j$ • x J 
i=i V i=l / 



(F4) 



i=i 

Performing the integral over x, we are left with the following integral over the tu/'s: 



(I x. m /• m I ra \ 

—J J J] exp f - £ ■ A- 1 • y M /4j 

/ i x. m „ m / m \ 

= [det (A/vr)]- 1/2 (-) J J] exp f - £ j (F5) 



The gaussian integration over the ui's leads to the final result 



[det (A/tt) det (B/47r)]- 1/2 (^)" (F6) 
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FIGURES 




u = cos 

FIG. 1. The tangent oscillates c — > b — > c — > b — > c. The curve is composed of two symmetric 
parts. 




FIG. 2. The family of curves ranges from the circle in the XY plane to the figure eight in the 
YZ plane. The perspective is slightly asymmetric to aid visualisation. The "fins" on these, and 
all other pictured curves, trace out the imbedded x and y axes, and thus depict the twisting of the 
rod. 
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writhe 



FIG. 3. The three types of behaviour of Ck for our writhing family of curves. The ratio A/C 
controls the plectonemic transition. 
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FIG. 4. The energy as a function of writhe for in the cases plotted in Figure ||. The energy is 
expressed in units of the combination 2A/1, where I is the circumference of the writhing loop, and 
A is the bending modulus. Also shown in this plot is the limiting value of the energy, at the figure 
eight configuration, which is displayed as a horizontal line. 
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FIG. 5. Typical potential, as given by Eq. ( A10| ), in the Schrodinger-like equation, (|A9|), obeyed 
by fluctuations about the deformed state. The quantity L is the total arclength of the loop. 
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FIG. 6. The translational mode, <&t(s)- The normalization of the mode is such that its derivative 
is unity at s = 0. The quantity L is the total arclength of the loop. 
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Stability determinant 

150 
100 
50 



■0.6 



■0.2 



-50 
-100 

FIG. 7. The determinant of the four-by-four matrix obtained by sandwiching the response 
kernel for the deformed loop between four of the five functions associated with the constraints 
on fluctuations about a solution of the Euler-Lagrange equations for a closed loop. Here the 
parameter a, b and c are equal to 1.14903, 0.145102 and -0.154898, respectively. The poles in the 
determinant are at the locations of the negative A's associated with unconstrained fluctuations 
about the extremal solution. 



Stability determinant 

10 




-1 -0.6 -0.2 

FIG. 8. A more detailed version of Fig |7|, in which the region between the pole at the less 
negative A and the origin is displayed, to highlight the fact that the determinant does not go 
through zero when A is negative. 
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FIG. 9. The determinant of the one-by-one matrix obtained by sandwiching the response ker- 
nel for the deformed loop between the function associated with the constraint of closure in the 
y-direction. The parameters a, b and c are equal to 1.14903, 0.145102 and -0.154898, respectively. 
The pole in the determinant is at the location of the negative A associated with an unconstrained 
fluctuation about the extremal function. 




FIG. 10. The deformation with fivefold symmetry. This represents one of the most "developed" 
member of the family of fivefold deformed loops. The perspective is off-center to aid the visualiza- 
tion. The actual curve is fivefold symmetric. The "fins" on this curve trace out the imbedded x 
and y axes, and thus depict the twisting of the rod. 
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FIG. 11. A trefoil solution to the Euler-Lagrange equations. This is a member of a family of 
solutions to the Euler Lagrange equations that knots in this particular way. It represents one of 
the more strongly "writhed" members of this particular class of solutions. The "fins" on this curve 
trace out the imbedded x and y axes, and thus depict the twisting of the rod. 



FIG. 12. A fivefold knotted solution to the Euler-Lagrange equations. This is a member of a 
family of solutions to the Euler Lagrange equations that knots in this particular way. It represents 
one of the more strongly "writhed" members of this particular class of solutions. The "fins" on 
this curve trace out the imbedded x and y axes, and thus depict the twisting of the rod. 
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FIG. 13. The function T(A), as denned by Eq. (|E6| ) for a characteristic supercoiling solution 
to the Euler-Lagrange equations for a closed loop. Note the zeros of this function at negative 
values of A. The function also passes throught zero at A = as the result of the existence of the 
translational mode. 
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FIG. 14. A more detailed plot of T(A) in the vicinity of A = 0. This function passes through 
zero with a positive slope, and there is, in addition, a zero at small, positive, A. 



39 



